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USE  OF  SURFACE  CHEMKIN  TO  MODEL  MULTIPHASE 
ATMOSPHERIC  CHEMISTRY:  APPLICATION  TO  NITROGEN 

TETROXIDE  SPILLS 

BRIAN  B.  BRADY  and  L.  ROBBIN  MARTIN 

Mechanics  and  Materials  Technology  Center,  The  Aerospace  Corporation,  P.O.  Box  92957,  Los  Angeles, 

CA  90009,  U.S.A. 

[First  received  17  March  1994  and  in  final  form  15  September  1994) 

Abstract— SURFACE  CHEMKIN  is  a  widely  available*  computer  program  developed  for  kinetic 
modeling  of  chemical  vapor  deposition.  We  show  that  it  may  be  adapted  for  kinetic  modeling  of  multiphase 
chemistry  in  the  atmosphere,  with  broad  capability  to  deal  with  complex  chemistry  and  physics.  It  can  deal 
with  multiple  phases  having  different  reaction  manifolds  in  each  phase,  it  deals  with  gas,  surface,  and  bulk 
reactions  and  mass  transfer  rates,  it  keeps  track  of  the  phase  equilibria  with  realistic  activities,  and  it  can 
operate  in  an  adiabatic  mode  to  include  the  effect  of  heat  release  on  the  system.  The  adapted  model  is 
applied  here  to  the  problem  of  a  nitrogen  tetroxide  spill  in  the  troposphere.  The  model  is  able  to  predict  the 
formation  of  a  nitric  acid /water  aerosol  and  to  follow  the  chemistry  taking  place  in  both  the  gas  and  liquid 
phases  as  the  spill  dilutes  in  the  surrounding  atmosphere.  The  model  predicts  that  in  such  a  spill,  most 
(70-90%)  of  the  nitrogen  oxides  released  are  converted  to  nitric  acid  over  a  wide  range  of  relative  humidity. 

Key  word  index:  Heterogeneous  chemistry,  kinetics,  dispersion,  toxic  plume. 


INTRODUCTION 

The  importance  of  heterogeneous  or  multiphase 
chemistry  in  the  atmosphere  has  only  recently  been 
recognized.  Chemistry  in  the  aqueous  phase  of  cloud 
drops  was  shown  to  be  important  in  the  atmospheric 
transformation  and  transport  of  sulfur  in  connection 
with  acid  rain  in  the  early  1980s  (Martin,  1984;  Op- 
penheimer,  1983;  Daum  et  ai,  1984).  Surface  chem¬ 
istry  on  polar  stratospheric  clouds  has  been  shown  to 
be  of  key  importance  in  explaining  the  Antarctic 
ozone  hole  (McElroy  and  Salawitch,  1989;  McElroy 
et  aU  1986;  Solomon  et  ai,  1986).  Surface  chemistry 
has  also  been  shown  to  be  important  in  the  NO^ 
chemistry  of  the  troposphere  (Heikes  and  Thompson, 
1983)  and  in  the  mid-latitude  stratosphere  (Brasseur 
et  a/.,  1990;  Granier  and  Brasseur,  1992). 

Multiphase  chemistry  poses  special  problems  in 
computer  kinetic  modeling  for  a  number  of  reasons. 
(1)  Chemical  reactions  on  or  within  the  second  phase 
will  in  general  have  different  rates  from  those  in  the 
gas  phase,  so  more  than  one  manifold  of  rate  con¬ 
stants  may  be  required.  (2)  Mass  transfer  to  or  within 
the  second  phase  may  be  rate  limiting,  rather  than 


*  Recently  available  under  a  no-cost  license  agreement 
from  Sandia  National  Laboratories.  Contact  Fran  Rupley, 
Computational  Mechanics  Division,  8245,  Sandia  National 
Laboratories,  Livermore,  CA  94551-0969,  U.S.A.,  for  current 
distribution  procedure. 


chemical  kinetics.  (3)  The  formation  of  a  condensed 
phase  usually  involves  a  discontinuity  in  the  activity 
of  a  species,  which  is  difficult  for  computer  programs 
to  handle.  (4)  The  formation  of  a  condensed  phase  will 
in  general  release  a  great  deal  of  heat.  (5)  Concen¬ 
trated  spills  of  certain  gases  such  as  HCl  cause  aerosol 
formation  because  they  reduce  the  partial  pressure  of 
water  of  the  liquid  phase,  thus  permitting  liquid  to 
exist  at  less  than  100%  humidity  (deliquescence).  In 
such  cases,  the  amount  of  liquid  phase  present  will 
be  strongly  coupled  to  the  progress  of  the  chemical 
reactions. 

An  atmospheric  release  of  nitrogen  tetroxide  in¬ 
volves  all  of  these  problems.  Upon  release,  the  follow¬ 
ing  sequence  of  events  will  take  place:  first,  there  is 
rapid  but  incomplete  dissociation  of  the  nitrogen  tet¬ 
roxide  (N2O4)  to  form  nitrogen  dioxide  (NO2).  Both 
the  dioxide  and  the  tetroxide  react  at  different  rates 
with  small  amounts  of  water  adsorbed  on  the  surface 
of  naturally  occurring  aerosols,  leading  to  the  forma¬ 
tion  of  nitric  acid,  which  is  deliquescent.  Water  vapor 
in  the  ambient  air  will  condense  in  the  nitric  acid, 
producing  a  liquid  phase  in  increasing  amounts.  This 
liquid  phase  is  a  more  favorable  medium  for  chemical 
reactions  between  the  oxides  of  nitrogen  and  water,  so 
the  entire  process  will  accelerate.  Thus,  a  cloud  will 
grow  by  this  “bootstrap’’  process.  As  all  of  this  hap¬ 
pens,  the  heat  released  by  condensation  and  by 
exothermic  reactions  will  increase  the  temperature  of 
the  air  parcel,  and  this  in  turn  may  limit  the  formation 
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of  liquid.  Modeling  this  complex  process  is  a  chal¬ 
lenge. 

Several  models  have  appeared  in  the  literature 
which  deal  with  various  aspects  of  the  multiphase 
problem.  A  fog  chemistry  model  by  Pilinis  et  al.  (1992) 
includes  multiphase  chemistry  and  equilibria  as  well 
as  aerosol  dynamics,  mass  transfer  rates  and  wet  and 
dry  deposition.  Thermal  effects  and  mixing  with  the 
surrounding  atmosphere  are  not  included,  so  the 
model  was  not  designed  to  deal  with  the  unusual 
situation  of  a  concentrated  nitrogen  oxide  spill,  in 
which  the  effects  of  heat  release  and  dilution  cannot 
be  ignored.  A  model  by  Moller  and  Mauersberger 
(1992)  deals  with  both  liquid  and  gas-phase  chemistry, 
but  is  constrained  to  a  fixed  droplet  size.  Spill  models 
by  Balentine  and  Eltgroth  (1985),  Spicer  and  Havens 
(1987),  Spicer  et  ai  (1986)  and  Lyons  et  al  (1990)  deal 
with  the  physics  but  not  the  chemistry  of  the  dispers¬ 
ing  material.  A  model  by  Rodean  (1989, 1991)  and  the 
ADAM  model  (Raj  and  Morris,  1987;  Mullet  and  Raj, 
1990)  deal  with  the  formation  of  nitric  acid  from 
a  nitrogen  tetroxide  spill,  but  these  models  are  based 
on  thermodynamics  and  not  kinetics.  A  model  by  Cho 
et  al  (1992)  includes  detailed  multiphase  chemistry 
and  physics  as  well  as  thermal  effects,  but  is  designed 
to  look  primarily  at  the  microphysics  of  individual 
particles  or  droplets. 

In  this  paper,  we  show  that  a  recently  developed 
computer  chemical  kinetics  program  can  be  modified 
in  a  straightforward  way  to  deal  with  the  problem  of 
a  concentrated  spill.  The  program  is  SURFACE 
CHEMKIN,  developed  by  Robert  Kee  of  Sandia, 
Livermore,  and  available  from  the  laboratory  via  li¬ 
censing  agreements  (Coltrin  et  al.y  1991;  Kee  et  al, 
1991a,  b).  The  program  was  developed  to  model 
chemical  vapor  deposition,  but  lends  itself  to  complex 
atmospheric  chemistry  problems.  It  is  capable  of  deal¬ 
ing  with  multiple  phases,  with  different  reaction  mani¬ 
folds  in  each  phase;  it  deals  with  gas,  surface,  and  bulk 
reactions  and  mass  transfer  rates;  it  keeps  track  of  the 
phase  equilibria  with  realistic  activities;  and  it  can 
operate  in  an  adiabatic  mode  to  include  the  effect  of 
heat  release  on  the  system.  The  effects  of  dilution  of 
a  spill  with  time  can  also  be  included. 

We  have  applied  this  program  to  the  problem  of 
a  large  spill  of  nitrogen  tetroxide  in  the  atmosphere. 
This  problem  involves  all  of  the  complexities  de¬ 
scribed  above,  because  of  the  unique  properties  of  the 
nitrogen  oxide- water  system.  Nitrogen  oxide  spills 
are  a  potential  hazard  from  the  storage  and  handling 
of  rocket  oxidizer  (Watje,  1978;  Haas  and  Prince, 
1984),  and  experiments  have  been  conducted  on  the 
release  of  nitrogen  oxides  into  the  troposphere 
(McRae,  1985;  McRae  et  aL  1984,  1987). 


MODEL  DESCRIPTION 

The  computer  model  is  based  on  the  SURFACE  CHEM¬ 
KIN  code  developed  by  R.  J.  Kee  of  Sandia,  Livermore.  This 


package  contains  subroutines  for  calculating  thermodyn¬ 
amic  properties,  equilibrium  constants,  rates,  and  unit  con¬ 
versions,  for  parsing  character  data  such  as  species  and 
element  names,  and  reaction  mechanisms,  and  for  solving  the 
set  of  differential  equations.  SURFACE  CHEMKIN  is  a  set 
of  software  tools  and  a  subroutine  library  designed  to  help 
create  a  program  to  solve  a  specific  chemical  kinetics  prob¬ 
lem.  The  package  is  designed  to  solve  a  group  of  differential 
equations  subject  to  a  set  of  constraints;  various  properties 
of  the  system  may  be  extracted  along  the  way.  One  needs  to 
write  the  program  which  describes  the  chemical  system  of 
interest,  and  calls  the  appropriate  subroutines  to  do  the 
calculations.  Multiple  surface  and  bulk  phases  can  be  in¬ 
cluded.  Surface  rates  are  scaled  by  the  ratio  of  surface  area  to 
gas-phase  volume.  The  package  is  designed  to  model  chem¬ 
ical  vapor  deposition.  We  wrote  a  driver  program  incorpor¬ 
ating  some  changes  in  the  calculation  method,  added  some 
new  subroutines,  and  modified  some  of  the  original  routines. 
The  resulting  program  is  well  adapted  to  atmospheric  chem¬ 
istry.  Our  alterations  are  described  in  more  detail  below. 

The  plume  model  is  adiabatic.  Chemkin  calculates  the 
enthalpy  of  each  species  in  the  model.  As  species  are  created 
or  destroyed  through  reaction  they  add  or  subtract  energy 
from  the  plume;  mixing  with  the  air  also  effects  the  energy 
balance.  This  net  change  in  energy  is  divided  by  the  heat 
capacity  of  the  plume  to  obtain  the  temperature  change  for 
each  integration  step.  This  energy  conservation  equation 
is  included  as  a  governing  equation,  or  constraint,  for 
SURFACE  CHEMKIN. 

The  governing  equation  for  temperature,  alluded  to  above, 
is  based  on  conservation  of  energy,  and  is  included  here  for 
reference: 


dT 

dt 


^  r  dtf.  dn( 
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The  preceding  equation  describes  the  change  in  the  plume 
temperature  with  time.  The  numerator  in  this  equation  is  the 
total  enthalpy  change  for  the  plume  in  time  interval  dt.  The 
summation  is  over  species  i,  where  N  is  the  total  number  of 
species  considered;  here  N  is  25.  H"  is  the  enthalpy  in  mass 
units,  m,  is  the  molecular  weight  of  species  z,  and  «?**•*'  is  the 
number  of  moles  produced  by  gas,  g,  surface,  s,  or  liquid 
phase,  reactions  per  unit  volume  (area  for  phase  =  s);  the 
product  //J"  m,  dn^dt  represents  the  heat  gained  or  lost  to  the 
plume  due  to  the  production  (or  consumption)  of  a  particu¬ 
lar  chemical  species  in  gas-phase  reactions  in  the  time  inter¬ 
val  dt.  The  rate  of  species  production  from  surface  reactions, 
dn]/dt  is  normalized  by  4^,  the  ratio  of  surface  area  to 
gas-phase  volume.  Fr  is  the  ratio  of  liquid-phase  volume  to 
gas-phase  volume,  which  normalizes  the  rate  of  species  pro¬ 
duction  from  bulk  phase  reactions,  dnf.Jdt.  T  is  the  temper¬ 
ature,  <Tp  is  the  cross-sectional  area  of  the  plume,  and  is 
heat  capacity.  The  term  {T*"  —  T)do^lu^dt  CJ*"  represents 
the  heat  required  to  bring  ambient  air,  which  is  mixing  into 
the  plume,  to  the  same  temperature  as  the  plume  gases.  The 
superscript  air  refers  to  quantities  in  the  undisturbed  air 
outside  the  plume.  The  denominator  in  the  above  equation 
represents  the  specific  heat  of  the  plume  and  is  made  up  of 
the  specific  heat  of  the  gas  phase,  C*p,  and  condensed  phase, 
portions,  where  p  is  the  gas-phase  density,  and  the 
total  mass  in  the  liquid  phase  per  unit  volume  of  gas  is  m^. 
The  heat  capacity  of  the  surface  phase  is  negligible  and  is 
ignored. 

The  remaining  constraints  are  based  on  conservation  of 
mass  and  are  broken  down  into  separate  equations  for  each 
species  in  each  phase.  The  governing  equation  for  gas-phase 
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species  is 

d  Yi  m,  f  rnf  fV 

—  _  —  ^  ^ 
dt 


+  {Yf 


The  symbol  Y  is  used  for  the  gas-phase  mass  fraction.  The 
first  term  in  this  equation  is  the  molar  changes  due  to 
reaction  converted  to  a  mass  fraction,  and  the  second  term  is 
the  effect  of  dilution  on  plume  species.  This  equation  has 
a  similar  structure  to  the  numerator  of  the  first  governing 
equation. 

The  equation  for  surface  species  in  terms  of  site  fraction, 
Xj,  is  similar.  There  is  also  a  correction  for  coverage,  where 
the  density  of  surface  sites  in  a  given  surface  phase;  is  given 
by  Fj  and  the  number  of  sites  occupied  by  a  given  species  is 
CT,.  Large  molecules  could  occupy  several  sites,  but  all  the 
species  in  this  model  are  small  molecules  with  a,  =  1.  Only 
one  surface  phase  is  used  in  our  model.  There  is  an  equation 
like  the  following  one  for  each  surface  species: 


dt  \  dt  J  Fj 


The  equations  for  liquid  species  are  also  similar,  but  given  in 
mass,  Zj.  A  term  for  dry  deposition  is  added.  Multiple  bulk 
phases  are  possible,  but  only  one  is  used  here.  An  equation 
like  the  following  is  generated  for  each  species  in  the  liquid 
phase. 


dt 


+  {zf 


where  the  final  term  is  an  approximate  correction  for  dry 
deposition  described  below.  This  complete  set  of  constraints 
is  used  in  solving  the  differential  equations  from  the  kinetics 
for  each  time  step  in  the  simulation. 

We  have  made  several  changes  to  the  Sandia  package.  To 
model  an  aerosol  we  keep  track  of  the  total  condensed 
volume  and  let  the  surface  area  vary  with  the  volume  and 
droplet  size.  We  created  a  separate  set  of  pure  bulk  phase 
(non-surface)  reactions  whose  rates  are  scaled  by  the  ratio  of 
condensed  phase  volume  to  gas  volume.  We  also  needed  to 
add  the  effect  of  atmospheric  dispersion.  For  purposes  of  this 
demonstration,  we  use  a  simple  Gaussian  plume  model,  and 
we  consider  only  the  centerline  concentration: 


r  \nayO./ 


where  R  is  the  evaporation  rate  of  N2O4,  r  is  the  wind  speed, 
and  Sy  and  s,  are  the  plume  spreading  parameters  for  the 
horizontal  and  vertical  directions.  For  simplicity,  the  spread¬ 
ing  parameters,  (Jy  and  cr.,  are  taken  to  increase  linearly  with 
time:  a,  =  <t,(0)  +  ax;  (7,(0)  +  hx;  for  this  case,  a,  =  1 
+  0.07x  (m)  and  =  2  +  O.lOx  (m),  which  corresponds 
roughly  to  the  Pasquil!  stability  class  “C’;  this  essentially 
means  that  in  the  absence  of  reactions,  the  concentration 
drops  with  the  square  of  time.  The  time  derivative  of  this 
equation  is  included  in  the  model.  More  sophisticated  dis¬ 
persion  models  are  available,  but  the  purpose  of  this  test  was 
primarily  to  see  the  effect  of  the  aerosol  chemistry. 

A  second  set  of  initial  conditions  is  used  for  the  atmo¬ 
spheric  temperature  and  composition;  background  atmo¬ 
spheric  gases  in  the  plume  increase  at  the  same  rate  the 
plume  gases  decrease.  The  total  pressure  is  conserved.  The 
program  was  run  at  this  point  with  a  helium  tracer  added  to 
the  initial  plume.  The  helium  mole  fraction  is  always  less 
than  !%  and  thus  does  not  appreciably  alter  the  density  of 
the  plume.  The  helium  does  not  participate  in  reactions, 
therefore  its  concentration  decreases  strictly  with  the  square 
of  time.  This  allowed  us  to  test  separately  the  effects  of 
reaction  and  dilution.  Conservation  of  mass  requires  that  the 


sum  of  all  moles  of  nitrogen  oxides  in  all  phases  divided  by 
the  moles  of  helium  be  a  constant: 

Z  ”nO  +  ”nO:  +  ”hN02  +  ”hNO, 

P  =  8.s.aq  _ _ _ _ _ ^ 


Two  more  changes  have  been  added  to  the  program. 
Approximate  corrections  for  dry  deposition  of  aerosol  and 
for  photolysis  rates  are  included.  For  dry  deposition,  we 
make  the  following  correction  to  the  concentrations,  based 
on  the  source  depletion  model  (Overcamp,  1976;  Horst, 
1977;  Corbett,  1981): 

C  =  Co[a,(0)/a,(r)]^'‘‘^'^'-“ 

where  Tj  is  the  deposition  velocity,  a  is  the  coefficient  of  x  in 
the  expression  for  as  shown  above,  and  u  is  the  wind 
speed.  A  nominal  dry  deposition  velocity  of  1  cms  ‘ 
(O.OI  Ms'M  is  chosen  for  NO2  (McRae  and  Russell,  1984) 
and  for  HNO3  (Huebert  and  Robert,  1985)  and  for  the  liquid 
droplets,  since  the  settling  velocity  for  the  model  calculated 
droplets  is  also  roughly  1  cms'L  Photolysis  rates  are  en¬ 
tered  as  the  average  daytime  rate.  Since  the  plume  is  at  first 
optically  thick,  the  overall  photolysis  rate  was  approximated 
by  multiplying  by  ///q,  where 

///o= 

Little  difference  in  the  final  results  was  seen  with  ///o  equal 
to  0  or  1,  so  this  approximation  of  reduced  but  uniform 
photolysis  is  not  critical:  The  adsorption  constant  for  NO2, 
a,  was  taken  as  15cm“^atm■^  P  is  the  pressure  of  NO2. 

The  bulk-phase  activity  coefficients  were  set  equal  to  the 
liquid  mole  fraction  with  two  exceptions.  Water  and  nitric 
acid  form  an  azeotrope;  the  real  activity  coefficients  for  these 
species  were  fitted  to  a  quadratic  equation  in  the  mole  fraction 
and  this  equation  was  used  in  the  program.  Thermodynamic 
data  for  the  bulk  species  not  available  in  SURFACE 
CHEMKIN  were  taken  from  the  JANAF  Thermochemical 
Tables  (Chase,  1985)  and  The  International  Critical  Tables 
(National  Research  Council,  1928);  surface  species  are  as¬ 
sumed  to  have  the  same  properties  as  bulk  phase  species  for 
this  simulation.  The  temperature  dependence  of  the  gas- 
phase  rates  is  included. 


MODEL  TESTS 

Our  first  lest  of  the  model  involved  simple  conden¬ 
sation.  We  created  a  hypothetical  gas-phase  reaction 
that  ran  at  room  temperature  and  had  water  vapor  as 
a  product.  The  model  correctly  determined  that  the 
gas-phase  w'ater  concentration  increased  until  it 
reached  the  equilibrium  vapor  pressure,  at  which 
point  liquid  water  began  to  form.  A  temporary  over¬ 
shoot  in  the  water  vapor  pressure  could  occur  if  the 
initial  reaction  rate  was  much  faster  than  the  rate  of 
condensation. 

A  second  test  of  the  program  was  for  the  case  of 
constant  liquid  to  gas  volume  ratio  and  constant 
droplet  size  because  we  were  able  to  solve  this  case 
with  a  hand  calculation.  We  let  dissolved  NOj  react 
with  liquid  water, 

2N02(aq)  +  H20(1)— ►HN03(aq)  +  HN02(aq) 

but  set  the  rate  for  evaporation  or  condensation  of  the 
water  to  zero.  The  half-life  of  the  NO2  is  calculated  in 
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this  case  as 

1 

where 

kV^H^ 

X'  = - 5 - . 

(1  +  f^RH) 

The  quantity,  Kr,  is  the  ratio  of  the  liquid-phase 
volume  to  gas-phase  volume,  H  is  the  Henry’s  law 
constant  for  NO2  in  water,  C®  is  the  initial  gas-phase 
concentration  of  NO2,  and  k  is  the  reaction  rate 
constant.  We  also  checked  the  mass  balance  to  insure 
that  there  was  no  net  loss  or  creation  of  material.  For 
atmospheric  reactions  the  pressure  is  held  constant. 

The  model  was  constructed  to  simulate  the  effects 
of  an  accidental  release  or  spill  of  liquid  nitrogen 
tetroxide,  N2O4,  or  NTO,  which  is  used  as  a  rocket 
propellant  in  large  launch  vehicles.  The  initial  plume 
from  the  spill  is  assumed  to  be  99%  gaseous  NTO  and 
1  %  helium.  The  helium  is  for  diagnostic  purposes  as 
explained  earlier.  The  effluent  expands  and  mixes  with 
air  in  accordance  with  a  Gaussian  model.  The  wind 
speed  is  10  km  h“  ^  and  the  dispersion  is  an  intermedi¬ 
ate  case.  The  background  air  is  made  up  of  nitrogen 
and  oxygen  in  the  usual  abundances,  20  ppb  nitrogen 
dioxide,  and  water  vapor  varying  with  the  humidity. 
The  background  air  is  also  assumed  to  contain  a  mar¬ 
ine  aerosol,  0.2  ppbv,  half  a  micron  in  diameter,  and 
consisting  of  water  and  sodium  chloride.  Such  an 
aerosol  is  typical  at  coastal  sites,  such  as  Cape 
Canaveral  and  Vandenberg  Air  Force  Base.  The  pro¬ 
portions  of  water  and  salt  in  the  aerosol  vary  with  the 
humidity.  For  simplicity,  aerosol  particles  are  also 
assumed  to  be  present  in  the  initial  NTO  plume.  The 


simulation  was  run  at  1  atm  pressure  and  293  K  for 
1  h. 

The  reaction  mechanism  includes  12  gas-phase  re¬ 
actions,  1 1  surface  reactions,  and  two  bulk-phase  re¬ 
actions  and  is  illustrated  schematically  in  Fig.  1.  In 
this  figure  the  boxes  represent  species  within  a  speci¬ 
fied  phase.  Arrows  connecting  the  boxes  represent 
chemical  reactions  or  mass  transport  paths  connect¬ 
ing  the  species,  with  the  principle  reactants  indicated 
beside  the  arrows.  A  complete  list  of  the  reactions  and 
rate  constants  are  listed  in  Table  1.  Mass  transport 
between  phases  in  handled  via  surface  reactions.  Mass 
transport  is  assumed  not  to  be  rate  limiting  except  for 
the  case  of  reaction  of  N2O4  with  the  droplet  surface, 
as  discussed  below.  Using  a  formula  for  the  character¬ 
istic  time  for  macroscopic  mass  transport  (interstitial 
gas-to-droplet  diffusion)  developed  by  Schwartz 
(1984)  (Levine  and  Schwartz,  1982): 

where  Dg  is  gas-phase  diffusivity,  L  is  liquid  volume 
fraction,  and  D  is  the  droplet  diameter.  For  our  cal¬ 
culated  cloud  properties,  we  get  a  characteristic  time 
of  about  0.1  s.  Since  the  model  results  without  mass 
transport  limitation  show  that  most  of  the  aerosol 
chemical  reactions  take  place  on  a  scale  of  about  10  s, 
macroscopic  mass  transport  should  not  be  a  limita¬ 
tion. 

The  rates  of  the  gas-phase  reactions  are  taken  from 
the  NIST  Database  (1992),  and  the  aqueous  phase 
rates  from  Park  and  Lee  (1988).  The  decomposition  of 
nitrous  acid  has  been  proposed  to  be  a  surface  reac¬ 
tion  by  Kobayashi  et  al  (1976),  and  the  program 
could  easily  incorporate  this  kind  of  reaction  if  the 
surface  reaction  rate  could  be  quantified.  For  the 


Source 


Fig.  1.  Schematic  diagram  of  the  main  chemistry.  The  boxes  represent  species 
within  a  specified  phase.  Arrows  connecting  the  boxes  represent  chemical  reac¬ 
tions  or  mass  transport  paths  connecting  the  species,  with  the  principle  reactants 
indicated  beside  the  arrows.  A  complete  list  of  reactions  is  given  in  Table  1. 
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Table  1.  Reaction  mechanism 


Reactions 

Rates 

Gas  phase 

O3  +  NOj^riNOj  +  O2 

1.14X  10” 

0  +  NOj;=:02  +  NO 

3.98x 

0  +  N02  +  M:^N03-f  M 

5.32  X  10” 

0  +  0, +  M;:i03  +  M 

3.35xlO’^e^«"/^ 

OH  +  NO2  +  M  ^HNO,  +  M 

3.04  X  10'  e””'’' 

N0-hN03;z±2N02 

6.32  X  10” 

NO  +  03^N02  +  O2 

1.47  X  10” 

NO2  +  NO3  +  M  ^N205  +  M 

2.86  X  10” 

N2O5  — NO2  +  NOj 

1.59  X  10”  e" 

N2O4  +  M^2N02  +  M 

1.46x  10”e'”“‘'’^ 

2H0N0;=±N02  +  no  +  H2O 

1.00  X  10”  6'“”^'’' 

2N0  +  02;=^2N02 

1.45x10” 

Surface 

N204(g)  -1-  H20(surf)^  -f  HONOjjgrf) 

0.001  <  <  1 

NaDsjg)  -|-  H20(5urf  Hl^D3(aq)  +  HN03(suff) 

v  =  1 

Bulk 

2N02(aq)  +  H20(/)^^  HN03(aq>  +  HONO,aq) 

8.4  X  10’ 

HONO,.„  +  H0N0,.„^N02,.„  +  NO,,,  +  H2O,.,, 

1.3  X  10“ 

Mass  transport 

+  H2G(surf>  ^  HaDjaurf)  ^20(1, 

}'=  1 

H2D(g)  +  NaCl(5urf)^=^  H20(surfj  NaCI(5q) 

=  1 

N02(g)  +  H20jsurr)^ND2(aq>  +  H20,5urf) 

r  =  1 

NO(g)  NO(aq)  H2D(syrf) 

V  =  1 

HNO3,,,  +  H20,.„rf,^=tHN03„,rf,  +  H2O,,, 

y  =  1 

MONO,,,  +  H20|,„,„;:=H0N0,.„,„  +  HjO,,, 

y  =  I 

HN03„„,r,  +  H20,„;^H20,.„.„  +  HNO3,.,, 

l.Ox  10” 

HON6,.„,r,  +  H20,„;=H20„„„  +  HONO,.„ 

1.0  X  10” 

NaCI(5yrf)  -f-  H20())  NaCljaqj  +  H20^surf) 

l.Ox  10” 

Note:  Gas-phase  reaction  rates  are  in  units  of  cm^  mol"’  s" ’  for  two-body  reactions 
and  cm^  mol  “  ^  s ”  ‘  for  three-body  reactions.  The  liquid-phase  rates  are  in  cm^  mol  ‘  ’  s  "  ’ 
and  must  be  multiplied  by  (0.0555)^  to  be  used  correctly  by  the  the  program.  All  gas-phase 
rates  are  from  NIST  (1992)  and  the  liquid-phase  rates  from  Park  and  Lee  (1988).  The 
gas-surface  reactions  and  the  surface-bulk  reactions  are  taken  not  to  be  rate  limiting 
except  for  the  N2O4  -h  water  reaction,  as  described  in  the  text.  T  is  the  temperature,  y  is 
the  sticking  coefficient.  The  species  subscripts  refer  to  the  phase:  (g),  gas,  (surf),  surface, 
(aq),  aqueous. 


present,  however,  we  will  include  it  as  a  bulk  aqueous 
phase  reaction  with  the  rate  given  by  Park  and  Lee 
(1988).  The  rates  of  gas-surface  reactions  are  put  in  as 
sticking  coefficients,  with  all  of  these  equal  to  1  except 
that  for  the  N2O4  +  water  reaction,  which  was  varied 
from  1  to  10"^,  and  0  to  “turn  ofT“  the  heterogeneous 
chemistry.  The  surface-bulk  reactions  are  put  in  as 
diffusion  limited.  Fifteen  gas-phase  species  are  in¬ 
cluded,  as  well  as  four  surface  species  and  six  bulk 
species.  The  program  tracks  the  concentration  of  each 
species  in  the  plume,  as  well  as  the  temperature,  den¬ 
sity,  particle  radius,  surface  to  volume  ratio,  liquid 
volume  fraction,  and  time.  For  our  conditions,  the 
program  runs  in  less  than  15  min  on  an  IBM  RISC 
System/6000  workstation. 


RESULTS  AND  DISCUSSION 

The  results  of  a  typical  run  are  shown  in  Fig.  2.  For 
this  run,  the  conditions  were  as  follows;  wind 


speed  =  lOkmh’^  source  strength  =  2(X)  kg  of 
N2O4  per  minute,  with  1  %  helium  included  as  a  refer¬ 
ence;  Gaussian  plume  dispersion  corresponds  rough¬ 
ly  to  Pasquill  stability  class  “C'’,  pressure  =  1  atm; 
temperature  =  293  K;  background  NO 2  =  20  ppb; 
background  aerosol  =  0.2  ppbv  of  sodium  chloride 
aerosol,  0.5  /tm  in  diameter;  water  vapor  in  accord¬ 
ance  with  the  relative  humidity,  which  is  95%  in 
Fig.  2;  and  sticking  coefficient  of  N2O4  =  1. 

It  can  be  seen  in  the  figure  that  the  volume  of  aerosol 
reaches  a  maximum  at  about  5  s  downstream  of  the 
release,  and  most  of  the  chemical  reactions  take  place 
between  2  and  30  s  downstream.  Depletion  of  the 
atmospheric  water  vapor  and  a  temperature  rise  may 
be  seen  for  a  brief  period  around  5  s.  After  this  time, 
mixing  of  the  ambient  air  restores  the  humidity  and 
evaporates  the  aerosol.  Note  that  the  acid  concentra¬ 
tion  within  the  droplets  falls  as  the  plume  dilutes  even 
though  the  droplets  are  evaporating.  We  also  note 
that  in  spite  of  the  temperature  rise,  at  no  time  does 
the  cloud  become  lighter  than  the  surrounding  air. 
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Adiabatic  Parcel  Evolution 

Fig.  2.  Time  evolution  of  a  parcel  of  the  plume  with  average  properties.  The  time  base  is  logarithmic.  Properties  shown  are: 
aerosol  radius;  liquid: gas  volume  ratio  (logarithmic);  temperature  (linear);  liquid-phase  (i.e.  within  aerosol  droplet) 
concentration  of  water,  nitric  acid,  and  nitrous  acid;  gas  phase  NO2,  N2O4  and  He;  gas-phase  water  (all  logarithmic 

ordinates). 


In  order  to  observe  the  effects  of  humidity,  we  ran  times  the  relative  humidity.  The  peak  temperature 

the  same  simulation  at  95,  80,  65,  50,  30,  20,  and  10%  and  liquid  volume,  and  the  final  gas  phase  NO2 

initial  relative  humidity.  The  gas-phase  water  concen-  pressure  and  liquid  mole  fraction  of  water  are  plot- 

tration  is  the  vapor  pressure  of  pure  water  at  293  K  ted  in  Fig.  3  as  a  function  of  humidity.  The  peak 


Humidity  dependence  of  Model 

Fig.  3.  Summary  of  runs  with  differing  initial  relative  humidity.  Shown  are  fraction  of  NO2  removed  by  aerosol  chemistry; 
peak  temp>erature  rise;  water  content  of  aerosol  after  1  h;  and  peak  liquid; gas  volume  ratio.  Sticking  coefficient  of  N2O4  is 

unity. 


Dependence  on  Sticking  Coefficient 

Fig.  4.  Summary  of  runs  with  differing  values  of  the  sticking  coefficient  of  N2O4  on  liquid  water.  Same  properties  are 

shown  as  in  Fig.  3.  Initial  relative  humidity  is  95%. 
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temperature  and  liquid  volume  increase  linearly  with  a  relative  humidity  of  95%  (after  correcting  for  the 
humidity,  whereas  the  final  N02(g)  and  water  concen-  background  NO2). 

trations  asymptotically  approach  a  limit.  At  10%  We  also  investigated  the  effect  of  the  surface  re¬ 
humidity  we  already  see  80%  conversion  of  the  nitro-  action  of  N2O4  by  varying  the  sticking  coefficient 
gen  oxides  to  nitric  acid,  and  this  rises  to  98%  at  from  one  to  zero.  The  results  are  plotted  on  a  semi-log 


Isothermal  Parcel  Evolution 

Fig.  5.  The  results  of  an  isothermal  run.  An  isothermal  situation  is  not  realistic  for  a  dense  spill  such  as  this,  but  it  illustrates 
the  effects  of  heat  on  the  plume  evolution.  See  text  for  discussion. 
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Fig.  6.  This  log-log  plot  shows  the  NO2  concentration  as  a  function  of  time  with  and  without  aerosol 
chemistry.  The  surface  reaction  of  N2O4  is  turned  oflf  by  setting  the  sticking  coefficient  to  zero. 
Aqueous  reactions  are  turned  off  by  setting  the  rates  to  zero. 


scale  in  Fig.  4.  The  water  content  of  the  aerosol  is 
insensitive  to  the  sticking  coefficient.  The  amount  of 
conversion  to  nitric  acid  is  very  sensitive  to  this  para¬ 
meter,  but  a  comparatively  large  number  (10"^  or 
higher)  seems  likely  in  view  of  earlier  studies  (England 
and  Corcoran,  1974).  The  other  relationships  are 
more  complex. 

The  effect  of  adiabatic  modeling  can  be  seen  in 
Fig.  5,  which  presents  the  results  for  an  isothermal 
situation  for  comparison.  An  isothermal  situation  is 
not  physically  realistic  for  a  concentrated  spill,  but  the 
comparison  gives  insight  into  the  effect  of  heat.  Inter¬ 
estingly,  there  is  not  a  great  difference  in  the  amount 
of  final  conversion  of  nitrogen  oxides  to  nitric  acid. 
These  results  were  initially  surprising  to  us,  since  we 
anticipated  that  the  temperature  increase  due  to  reac¬ 
tion  and  condensation  would  greatly  reduce  the 
amount  of  aerosol  chemistry.  The  comparatively 
small  effect  on  the  outcome  from  including  adiabatic 
effects  appears  to  be  due  to  the  initial  temperature 
decrease  resulting  from  the  evaporation  and  dissocia¬ 
tion  of  the  N2O4.  This  initial  decrease  in  temperature, 
which  is  ignored  in  the  isothermal  treatment,  causes 
less  dissociation  of  the  N2O4  to  NO2,  and  conse¬ 
quently  a  more  rapid  initial  reaction  of  the  N2O4  with 
the  aerosol.  Thus,  the  initial  fall  in  temperature  offsets 
the  subsequent  rise  in  temperature. 

Figure  6  compares  the  results  with  and  without 
aqueous  (bulk)  and  surface  chemistry,  i.e.  with  and 
without  the  effect  of  aerosol  chemistry.  According  to 
the  model,  inclusion  of  the  aerosol  chemistry  reduces 


the  amount  of  gas  phase  NO2/N2O4  by  about  a  fac¬ 
tor  of  ten.  (This  ratio  decreases  after  about  an  hour 
due  to  mixing  in  of  fresh  background  NO2.)  The 
balance  is  converted  to  nitric  acid  which  will  partly 
rain  out  as  liquid,  and  partly  evaporate  into  the  gas 
phase. 

SECONDARY  CHEMISTRY  IN  THE  PLUME 

In  a  separate  set  of  computer  runs,  we  included  an 
expanded  set  of  reactions  designed  to  mimic  the  essen¬ 
tial  chemistry  leading  to  the  formation  of  ozone  and 
OH  radicals.  The  reason  for  our  interest  in  this  chem¬ 
istry  is  that  if  there  were  a  substantia!  increase  in  the 
concentration  of  OH  radicals,  then  there  would  be  an 
increased  oxidation  rate  of  hydrocarbons  and  other 
carbon  containing  molecules,  which  would  lead  in 
turn  to  increased  formation  of  peroxyacetyl  nitrate 
(PAN).  PAN  is  an  important  nitrogen  containing  spe¬ 
cies  in  urban  air  pollution. 

The  reaction  set  added  to  the  basic  model  is 
as  follows  (k's  are  effective  second-order  rates,  cm^ 


molecule  's  '): 

NO2  +  hr  =  NO  +  0(’P), 

9x  10"^  s'* 

0(^P)  +  O2  =  O3, 

1.5x10'*“ 

O3  +  hv  =  0(*P)  +  O2, 

3.6xl0'‘‘s' 

=  0('D)  +  O2, 

1.6  X  10' 's' 

0('D)  +  M  =  0('P)  +  M, 

2.9  X  10'** 
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Fig.  7.  This  plot  depicts  the  effect  of  the  plume  on  OH  concentrations.  OH  concentration  is 
plotted  on  a  log  scale  vs  time  on  a  linear  scale  for  two  scenarios:  release  of  a  plume  of  N2O4,  and 
release  of  a  plume  with  the  same  composition  as  the  ambient  air.  The  OH  in  the  oxidizer  plume  is 
suppressed  by  two  orders  of  magnitude.  This  suggests  that  within  the  I  h  time  period  covered  by  the 
model,  hydrocarbon  oxidation,  and  therefore,  PAN  formation,  are  not  significant.  At  longer  times  the 
OH  should  recover,  leading  to  PAN  formation. 


0(*D)  +  H2O  -  20H,  2.2  X  10" 

OH  +  NO2  +  M  =  HNO3,  1.1x10" 

OH  +  NO  +  M  =  HNO2,  6.8  X  10"  " 

HNO2 -h /tv  =  OH  +  NO,  1.6xl0"^s"‘. 

The  results  of  adding  this  reaction  set  to  the  plume 
model  were  that  the  ambient  ozone  and  OH  concen^ 
trations  are  deeply  reduced.  Figure  7  shows  the  effect 
of  the  plume  on  OH  concentration  for  full  sunlight 
penetration,  which  gives  an  upper  limit  to  the  effect  of 
photochemistry.  Since  the  OH  concentration  is  re> 
duced  by  a  factor  of  over  100  during  the  modeling 
period  (1  h),  there  will  be  little  oxidation  of  hydro¬ 
carbons  and  therefore  little  production  of  PAN.  At 
longer  times  the  OH  concentration  will  recover,  and 
there  will  be  an  increase  in  PAN  due  to  the  increased 
NO,. 


CONCLUSIONS 

We  have  shown  that  the  computer  software  SUR¬ 
FACE  CHEMKIN  can  be  adapted  for  modeling  at¬ 
mospheric  aerosol  chemistry,  in  particular  the  situ¬ 
ation  of  a  nitrogen  tetroxide  spill.  The  model  here 


includes  the  effects  of  heat  release,  formation  and 
evaporation  of  an  aerosol,  and  gas,  liquid,  and  surface 
chemistry.  The  results  are  sensitive  to  a  reaction  not 
studied  in  the  literature,  namely,  the  surface  reaction 
of  the  gas-phase  molecule  N2O4  with  water  at  the 
aerosol  surface.  The  results  suggest  that  an  aerosol 
will  form  over  a  wide  range  of  humidities  (10-95%) 
and  that  the  resulting  chemistry  will  convert  roughly 
90%  of  the  nitrogen  oxides  to  nitric  acid. 

This  model  uses  a  simplistic  dispersion  model  but 
treats  the  multiphase  chemistry  rigorously.  If  a  more 
accurate  dispersion  model  were  included  it  could  pro¬ 
vide  an  even  better  picture  of  the  plume  dynamics,  but 
such  a  model  would  be  computationally  intensive.  In 
the  chemical  system  treated  here,  the  unique  aerosol 
chemistry  is  over  very  quickly.  This  suggests  that  the 
chemical  make-up  of  the  plume  at  100  s  calculated 
using  SURFACE  CHEMKIN  could  be  used  as  the 
source  term  in  a  more  precise  dispersion  model.  This 
would  give  a  more  accurate  concentration  down¬ 
stream  for  reasonable  computation  times. 
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